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ABSTRACT 

We have analyzed the first 3.75 years of data from TAOS, the Taiwanese American Occultation 
Survey. TAOS monitors bright stars to search for occultations by Kuiper Belt Objects (KBOs). 
This dataset comprises 5 x 10 5 star-hours of multi-telescope photometric data taken at 4 or 5 Hz. No 
events consistent with KBO occultations were found in this dat aset. We compute t he number of events 
expec te d for the Kuiper Belt formation and e voluti on models oflPan fc Saril (|2005l ) , iKenvon fc Bromlevl 
(|2004h . lBenavidez fc Campo Bagatinl (|2009t ). and lFraserl (|2009T >. A comparison with the upper limits 
we derive from our data constrains the parameter space of these models. This is the first detailed 
comparison of models of the KBO size distribution with data from an occultation survey. Our results 
suggest that the KBO population is comprised of objects with low internal strength and that planetary 
migration played a role in the shaping of the size distribution. 
Subject headings: Kuiper Belt, occultations, Solar System: formation 



1. INTRODUCTION 

The Kuiper Belt has been shaped by accretion and 
disruption processes throughout the history of the Solar 
System. With small orbital eccentricities, the relative 
velocities of the objects in the early Kuiper Belt were 
sufficiently low to allow accretion processes to form kilo- 
meter and much larger objects. Later, the velocity dis- 
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persion increased, possibly as the KBO population was 
stirred up by the gravitational effects of the larger plan- 
ets and planetoids. Only large objects were then able 
to continue growing through impacts, whereas collisions 
among smaller bodies resulted in disruption. The de- 
tails of these processes depend on the internal strength 
of the KBOs and on the orbital and dynamical evolution 
of the gas giant planets. The size distribution of KBOs, 
therefore, contains information on the internal structure 
and composition of the KBOs - and hence information 
on the location and epo ch in which they for med - and 
on planetary migration (jKenvon et al.l [20081 and refer- 
ences therein). Direct observatio ns have detected KBO s 
as faint as magnitude R ~ 28.2 (|Bernstein et al.l I2004T) . 
which corresponds to a diameter of about 27 km as- 
suming a 4% albedo. The large end side of the KBO 
size distribution can therefore be characterized through 
its brightness distribution. The latter is well described 
by a power law £(< R) = lQ"^--^), with an in- 
dex a = 0.6 and Rp = 23 (|Fraser fc Kavelaarsl [20091 
Fuentcs & Holman 2008) for objects brighter than about 
R = 25, or D ~ 100 km. This is the region of the size 
spect rum which reflects the early history of agglomera- 
tion. IKenvon fc Windhorst! (|2001[ ) pointed out that the 
intensity of the infrared Zodiacal Background sets limits 
on the extrapolation of a straight power law to smaller 
sizes. The relatively sh allow size distributio n of Jupiter 
Family Comets fJFCs. pTancredi et al.|[2006T) . which are 
believed to originate in the Kuiper Belt, and the cra - 
tering of Triton observed by Voyager 2 (jSternl fl996| ). 
all point to a flatter distribution for small KBOs 1 . In 
2004 evidence surfaced that a break in t he power law oc- 
curs a t a diameter larger than 10 km: iBernstein et al.1 
(2004) conducted deep Hubble Space Telescope obser- 
vations with the Advanced Camera for Surveys which 

1 The relationship between the cratering o f Triton and the 
Kuipe r Belt size distribution is questioned by Schcnk & Zahnld 
i2007h . 
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led to the discovery of only 3 new objects fainter than 
R = 26, about 4% of the number expected from a sin- 
gle power law distribution extrapolated to 10 km. While 
this work remains the state of the art for deep direct sur- 
veys of the Outer Solar System, recent campaigns have 
observed many more faint objects down to magnitude 
R = 27, which with the assumption o f a 4% albedo corre- 
sponds to about 40 km in diam e ter 2 dFraser fe Kavelaarsl 
20081 iFuentes fc Holmanll2008l iFuentes et all 120091 and 
Fraser fc Kavelaarsl 120091) . These recent data allowed 



TABLE 1 

Dataset parameters (3-telescope data) 



them to locate a break in the power law size distribu- 
tion at diameters 30 < D < 120 km. 

The region of the size spectrum between tens of kilo- 
meters and meters in diameter is particularly interest- 
ing as models predict here the occurrence of transi- 
tions between different regimes where the binding en- 
ergy of KBOs is dominated either by gravity or internal 
strength. These tr ansitions would leave a signature in the 
size distribution jP an & Sari 2005, Ken von fc Bromlevl 
12004 iBcnavid ez fc Campo Bagatinll2009l and references 
therein). Occultation surveys allow us to reach farther 
then the current limits of direct observations, and into 
this region of interest. These surveys monitor back- 
ground stars in order to detect the chance alignment of a 
KBO with a target star, which would generate a variation 
in the observed flux of the star. At distances in the Outer 
Solar System (tens to thousands of AU) the signature 
left in a lightcurve by the transits of D ~ 1 km objects is 
dominated by diffraction. This technique requires high 
frequency photometric time series as the time scale for an 
occultation by a n Outer Solar System Object is a frac- 
tion of a second faoques fc MoncuauetJl2000l iNihei et al.l 
120071 iBickerton et al. 2009). A few such surveys have 
been attempted in the past several ye ars and have 
recent ly started repo r ting r es ults: e.g.. Hloaues et al.l 
(I2006h. IChanget al.l (lj)07l). IBickerton et all (|2008f K 
iLiu et all ( 2008). Zhang et al. (2008) - h ereina fter Z08 - 
Bianco et alT(|2009D . and I Wang et all §009). None of 



these surveys have claimed detections in the Kuiper Belt; 
upper limits have thus been placed on the number den- 
si ty of KBOs in the sky. 

IBickerton et al.l J2008) set an upper limit to the sky 
density of KBOs of E N (D > 1 km) < 2.8 x 10 9 deg" 2 
using the 40 Hz dat a from their own su rvey as well as 
the 45 Hz data fromlRoques et al.l (120061) an d the X -ray 
data from lChang et all (|2007l ). iBianco et al.l (|2009D car- 
ried out a 30 Hz survey with Megacam at the MMT 
setting a more stringent limit of £at(-D > 1 km) < 
2.0 x 10 8 deg" 2 and a limit of T, N (D > 0.7 km) < 
4.8 x 10 8 deg" 2 . IWang et all (|2009l ) reported prelimi- 
nary analysis of videomode engineering data taken with 
the Pan-S TARRS system. 

Recently ISchlichting et al.l (|2009l ) reported the detec- 
tion of a candidate occultation event consistent with a 
D ~ 1 km KBO in the analysis of archival guiding data 
from HST, and an estimate of the sky density of KBOs 
of EjvCD > 0.5 km) = 2.1± 4 ; 8 x 10 7 deg~ 2 . 

The Taiwanese American Occultation Survey (TAOS) 
has been operating since 2005 with two, three, and now 

2 The magnitude of KBOs is converted into diameter by assum- 
ing a nominal 4% albed o throug hout the paper, note h owever that 
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Fraser & Kavclaars 2008 and Fraser & Kavclaars 2009 assumed an 
albedo of 6% in their work. 



a Multi-telescope measurements. 

four telescopes simultaneously taking stellar photometry 
at 5 Hz 3 . The analysis of the first two years of TAOS 
reported no detections (Z08) and an upper limit was 
derived to the slope of the small size end of the size 
spectrum. The TA OS system is described in detail in 
iLehner et al.l (2009). Using 50 cm aperture robotic tele- 
scopes in simultaneous observations and observing with 
the relatively low cadence compared to the aforemen- 
tioned occultation surveys, TAOS was designed to ad- 
dress the km-size region of the KBO size spectrum. We 
will show here that the marginal sensitivity to sub-km 
objects is more than compensated by the very large ex- 
posure of our star targets. Here we consider the first 3.75 
years of TAOS data, a significantly larger dataset than 
the one explored in Z08. With these data we are able to 
constrain Kuiper Belt formation and evolution models. 

In Section [2] we describe the new dataset. In Section 
[3] we briefly describe our detection algorithms, as well 
as our efficiency analysis. We also discuss our recovery 
efficiency and discuss the most productive strategies for 
TAOS and the other occultation surveys, and, in Sec- 
tion 13.41 we derive the effective coverage of our survey. 
In Section |4] we derive model-independent limits to the 
number of objects in the Kuiper Belt, and we compare 
our results with those of similar surveys. In Section[5]we 
briefly describe models for the formation and evolution 
of the Kuiper Belts and we then derive and discuss con- 
straints to these models. In Section [5] we compare our 
upper limits to the estimates on the number of KBOs 
set by dynamical simulations for JFC progenitor popu- 
lations. Finally we summarize and discuss our findings 
in Section 

2. 3.75 YEARS OF TAOS DATA 

TAOS is a dedicated survey that observes at a cadence 
of 5 Hz. The primary scientific goal of the survey is to 
estimate or set constraints on the number of KBOs in the 
region of the size spectrum that is currently too small to 
be observed directly: D < 10 km. 

Here we present an expanded analysis of three- 
telescope TAOS data 4 . These data consist of photo- 
metric measurements of target star fields collected syn- 
chronously with all three telescopes. The dataset ana- 
lyzed here was collected between January 2005 and Au- 
gust 2008. In a previous analysis of a subset of these 
data, Z08 reported an upper limit to the size distribu- 
tion of KBOs under the assumption of a single power 
law for small KBOs. If one models the size distribution 

3 A small subset of early data was collected at 4 Hz cadence, 
comprising about 5% of the data analyzed in this work. 

4 The fourth telescope, TAOS C, became operational in August 
2008. The results presented in this paper are based on analysis of 
all of the three-telescope data collected to this point. 
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Fig. 1. — Distribution of magnitudes for the TAOS target stars, (bin size 0.16 mag, left). SNR for the TAOS target stars, averaged over 
the duration of a run and over the three telescopes (bin size 0.73, right). A few targets at greater SNR and brighter magnitude, amounting 
to < 5% of the data, are not shown. 



for objects s maller than D = 28 k m, the smallest direct 
observation (jBernstein et al.|[2004h . as a single power law 
dN/dD oc D~ q , where N is the surface density of ob- 
jects, the slope of the distribution is limited to q < 4.6. 

Throughout the remainder of this paper, a data run 
refers to a set of data collected in an uninterrupted ob- 
servation of any field. For a single star in the field a set 
of three lightcurves belonging to one data run will be 
referred to as a lightcurve set, and each three-telescope 
measurement, at a single time point, will be referred to as 
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Fig. 2. — SNR versus TAOS instrumental magnitude AfxAOSi 
only a random sample of 1% of all stars is shown for clarity. 



a triplet. A star-hour refers to an hour of high-cadence, 
multi-telescope observations on a single target star. 

The data set described in this paper amounts to 
amount to 5.0 x 10 5 star-hours, while the data set used 
in Z08 comprises 1.5 x 10 5 star-hours. The details of this 
dataset, and of the dataset published in Z08, are sum- 
marized in Table [TJ Over 90% of our data is collected 
within 5° of the ecliptic plane in order to maximize the 
rate of occultations. 

TAOS uses the zipper-mode technique to read out the 
C CD cameras at high f requency. This method, described 
in iLehner et al.l (|2009| ). enables high speed observations 
across the 3D field of view of the TAOS telescopes, but 
it artificially increases the crowding of the field and the 
background. In zipper-mode readout each star in the 
field is represented in a subsection of the output image - 
which we call rowblock and which comprises 76 rows for 
our 5 Hz data - so that the field of view is entirely im- 
aged in each rowblock. Note that the images of different 
stars in a rowblock, however, do not necessarily belong 
to the same epoch. The zipper-mode readout boosts the 
sky background by a factor of 27 at a 5 Hz readout rate. 
This limits the sensitivity of TAOS to stars as faint as 
A^taos = 13.5, for which a signal-to-noise ratio (SNR) 
of ~ 7 can be achieved in a dark night. The magnitude 
and SNR distributions for the target stars in our survey 
are shown in Figure [TJ On the left panel, the x— axis 
is the TAOS instrumental magnitude MtaoSi which is 
defined by a regression on the USNO-B magnitudes to 
be similar to -Rusno- The correlation between instru- 
mental magnitude and SNR is shown in Figure [2] The 
scatter in the relationship between SNR and Mtaos is 
due to both changes in the sky background and in the 
weather conditions, and to different degrees of crowding 
in the fields. In Figure [3j left, we show the number of 
star-hours at different angles from opposition. The top 
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Fig. 3. — Left: distribution of angles from opposition for the TAOS targets. The top axis shows the relative velocity of a KBO at 43 AU, 
given the position of the field. The bin size is 10.5°. Right: distribution of ecliptic latitude for the TAOS target fields (center of the field 
is assumed), bin size ~ 2.5° 



scale indicates the velocity of a KBO at the center of a 
field at this elongation. We cover a large range of oppo- 
sition angles; our field selection algorithm favors ecliptic 
fields near zenith. Most angles are positive because the 
weather at the site tends to improve after midnight. The 
right panel of Figure [3] shows the distribution of ecliptic 
latitude of our data. The effects of the angle from op- 
position on our efficiency and event rate, as well as the 
efficiency as a function of magnitude and crowding are 
discussed further in Section [3.31 



3. ANALYSIS 

The first step in the analysis is the photometric reduc- 
tion of the zipper mode images in the data set. A custom 
aperture photometry package (|Zhang et al.ll2009l) is used 
to measure the brightness of each star at each epoch, and 
the resulting series of flux measurements are then assem- 
bled into lightcurves for subsequent analysis, which is 
described in the following subsections. In Section [37X1 wc 
describe our detection algorithm and the rejection of false 
positives. We then describe the efficiency tests: in Sec- 
tion [32] we describe in detail how we identify the angular 
size of our target stars to simulate occultations correctly. 
In Section [3.3l we show how we simulate and implant oc- 
cultation events in our lightcurves, and test the behavior 
of our efficiency as a function of various parameters rel- 
ative to the occultations and to the observing strategy. 
Finally, we can derive the effective coverage of our survey 



(Section 

3.1. Event Detection and False Positive Rejection 

The Fresnel scale is defined as F — (AA/2)^ where A is 
the wavel ength of observation , and A the distance to the 
occulter ([Rogues et all [19871: iBorn fc Wolfl [l980h . For 
optical observations at the distance of the Kuiper Belt 
(about 43 AU) the Fresnel scale is F « 1.4 km. Occul- 
tation events will therefore exhibit significant diffraction 



effects. Occultations are manifested in the lightcurve of 
an observed star as an alternation of bright and dark 
features, typically with an overall suppression of the 
flux. Theoretical occultation lightcurves are shown in 
Figure O The signature of an occultation by a KBO 
of sub-kilometer size has a duration of about 0.2 sec- 
ond at opposition, and about a second near quadrature. 
A typical KBO occultation is then expected to result 
in the suppression of the flux for one or a few consecu- 
tive points in a TAOS lightcurve. In order to ascertain 
the extra-terrestrial origin of a dip in a lightcurve TAOS 
observes simultaneously with multiple telescopes. This 
allows us to rule out, on the basis of simple parallax con- 
siderations, atmospheric scintillation phenomena which 
might mimic an occultation event and which could be a 
source of false positives in occultation surveys, as well 
as any non-atmospheric phenomena such as birds, air- 
planes, etc. 

In order to detect occultations we need to identify brief 
flux changes in a star simultaneously observed by all tele- 
scopes. The statistical significance of a simultaneous low 
point in our lightcurves can be assessed rigorously, and 
the probability of a low measurement being drawn out of 
pure noise decreases with the number of telescopes ob- 
serving the target, provided that the measurements for 
the telescopes are independent. The lightcurves are high- 
pass filtered to remove trends due to weather patterns 
and changes in atmospheric transparency. High-pass fil- 
tering the lightcurves preserves the information on time- 
scales relevant to occultation phenomena (one or a few 
points in a time series). The implementation of the filter 
is described in Z08. The filter produces a time series in 
which h(t), the measurement taken at time t, represents 
the deviation from the local mean of the lightcurve in 
units of local standard deviation. 

To detect events we rank-order the photometric mea- 
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Fig. 4.— Angular size distribution for a typical TAOS field (field 120, RA: 13°. 7, Dec: -10°. 7) derived from the 2MASS K — J colors 
(a). The curves show the theoretical behavior of the angular size for AO, F0, GO, K0 dwarf stars (thin lines) and GO, K0, K5 giants (thick 
lines). The size of the Fresnel scale at 43 AU is shown at 0.08 mas (dashed line). Best fit to the angular size distribution (b): x-axis is 
the USNO-B B — R color and y is the angular size derived using Equations [l] & [2] converting all stars to apparent magnitude R = 12. 
Implanted angular sizes, derived from USNO-B the B — R color for our simulation (c) and angular size of the star for which events are 
recovered (d); a random subset of 1% of our data is plotted in the bottom panels. 



surement in each of our lightcurves, from the lowest to 
the highest flux, independently for each telescope (la- 
beled A, B and D). The i-th point in a lightcurve will be 
associated to rank rj for telescope T. We then consider 
the rank triplets [rf, rf , rf). The probability distri- 
bution of the quantity Zi = — ln{r^ rf rf/Np}, with 
N p the number of points in the lightcurve set, can be de- 
termined combinatorially. Knowing this, under the null 
hypothesis that there is no event in the triplet i, we can 
compute the probability for a random variable Z aris- 
ing from this distribution P(Z > zi) = £. We set a 
threshold such that we expect fewer than 0.27 events in 
our dataset that are due to random fluctuations. For the 
dataset discussed in this paper we accept as events all 
data points that produce a rank product less likely than 
£ = 3.0 x 10~ n to be drawn from a random distribution. 
Note that events generated by large KBOs, or for ob- 
servations near quadrature, would affect more than one 
point in the lightcurve (Figure [5]) , and our rank-based 
search algorithm is most efficient when the dip in the 
lightcurve is isolated. Therefore, in addition to search- 
ing for single-point events, we also bin our lightcurves 
by 2, re-rank them and repeat the statistical tests de- 
scribed above. The probability of each data point is as- 
sessed for both unbinned and binned lightcurves. Each 
lightcurve is binned twice, with two different starting 
points. This increases the detectability of occultations 
by large KBOs and by KBOs transiting with low rela- 
tive velocit y. For a detailed dis cussion of our statistical 
analysis see lLehner et al.l (|2010l ). 

For a set of lightcurves of a given star, the ranks in the 
three telescopes should not be correlated for the statis- 
tical analysis described above to be valid. We have de- 
veloped a series of statistical tests to identify data runs 



where significant correlations (typically due to fast mov- 
ing cirrus clouds) are found in the lightcurves. In such 
data runs the ranks are not independently distributed 
and thus we can not accurately determine the statisti- 
cal significance of any candidate events. Any data run 
where the independence of the measurements after filter- 
ing cannot be rigorously established is removed from our 
dataset. A complete description of these tests is beyond 
th e scope of thi s pape r, but they are discussed in detail 
in lLehner et~aT1 (|2010D 

For the next step we relax considerably the ulti- 
mate significance requirement described above, and se- 
lect as provisional candidates those triplets that have 
£<1.0xl0~ 6 . Note that the significance £ refers to the 
probability that the point would be drawn from a ran- 
dom distribution, therefore the lower the value of £ the 
higher the statistical significance of the event. Nearly 
150,000 provisional candidates are found. We use all 
these measurements to identify and remove spurious re- 
gions of the lightcurves, and hence identify and reject 
false positive events which arise from sources other than 
random chance. The constraints described below allows 
us to recognize regions of lightcurves with atypical noise 
and contamination by transiting objects (satellites or me- 
teors, which turn out to be the major source of false 
alarms), and to identify high-frequency fluctuations in 
the raw data that are not removed by the high-pass fil- 
ter. These are the steps of our false positive rejection 
process: 

• Contiguity: Contiguous candidates within a 
lightcurve and candidates that are within three 
time-stamps of each other are removed. Only the 
one rank triplet that has the highest significance 
in a series of contiguous or proximate points is 
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Fig. 5. — Steps of the generation of a simulated occultation event. Top, the left and right panels both show the same point source 
lightcurve for a 3 km KBO at 43 AU occulting an F0V star. Second row: finite source lightcurves for the same occultation parameters 
for a V = 11 star, corresponding to an angular size of 0.015 mas, for a zero impact parameter (left), and at an impact parameter of 2 km 
(right). Row three: the lightcurves in row two are integrated over intervals of 105 ms for the occultation above at opposition (left) and at 
50° from opposition (right). The lightcurves are sampled at 5 Hz, with no time offset (left) and with a time offset of 50 ms (right). 



considered as a candidate. This removes double- 
counted events: events caused by large KBOs or 
KBOs moving at low relative velocity would affect 
more than one contiguous point. Furthermore this 
removes events that are double-counted because 
they appear significant in both the binned and un- 
binned lightcurves. This eliminates about 40% of 
the candidates. 

• Simultaneity: Candidates that appear in the 
lightcurves of more than one star simultaneously 
at the same time-stamp or within three rowblocks 
are considered to be false positives. We expect si- 
multaneous count drops in time-domain to be pri- 
marily due to inaccurate aperture positioning in 
the photometry. In the rowblock domain simulta- 
neous count drops might be due to inaccurate back- 
ground determination or non-occultation events al- 
tering the baseline of the lightcurve at, or around, 
the candidate event, or by fast moving cirrus clouds 
or other phenomena which induce high frequency 
fluctuations in multiple lightcurves which are not 
removed by the high-pass filter. This cut removes 
about 60% of the remaining candidates. 

• Number of telescopes: At this point we require all of 
our remaining candidates to have been observed by 
all three telescopes. Although we are only consid- 
ering three-telescope runs in our analysis, for some 
targets the lightcurve might not be extracted in the 
photometry phase for all telescopes. Small differ- 
ences in the field of view and in the field distortion 
might make one star target not visible to all tele- 
scopes if it is at the edge of the field or if the crowd- 
ing induced by the zipper-mode re adout caused 
overl ap of the target with other stars ()Lehner et al.l 
2009). About 35% of the remaining candidates are 



thus removed. Note that this cut cannot be applied 
earlier as simultaneous events might appear in only 
one 3-telescope lightcurve, but also in 2-telescope 
lightcurves, and we want to be able to recognize 
and remove these events. At this point there are 
still over 20,000 candidates left. Note that we do 
not count the star hours discarded by this cut in 
our total exposure of 5 x 10 5 star-hours. 

Significance threshold: We finally constrain £ such 
as to expect fewer than 0.27 false positives in our 
dataset. This constraint depends on the size of 
the dataset: for the 9 x 10 9 triplets remaining £ < 
3.0 x 10~ n allows < 0.27 false positives due to 
random noise. Only 228 candidates remain. 



The remaining candidates require visual inspection: 
first of the lightcurves, and for any remaining candi- 
dates, of the images. Most of the events are caused by 
the passing of bright objects, such as artificial satellites, 
meteorites or asteroids, that generate a variation in the 
background or baseline of the lightcurve responsible for 
causing artificially low counts in the neighborhood of the 
object. Many, but not all, of these false positives are 
removed by the simultaneity cut described above. Note 
that in our observing mode bright stars gene rate a bright 
strea k across the length of our images (see IZhang et al.l 
2009) , as flux is collected during the shutterless row shift- 
ing. In the presence of a bright object overlapping with 
a star-streak generated by the zipper-mode readout, the 
brightness of the streak is overestimated, thus too much 
flux is subtracted from the rowblock column causing an 
artificial flux drop in the star time series. In many in- 
stances the foreground object will also appear inside the 
star aperture artificially boosting its brightness. This 
flux drop will then be associated with a very high flux 
measurement following or preceding the event epoch, a 
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signature that allows us to remove these false positives by 
inspecting the lightcurve. We also inspect the centroid 
position of the aperture. If the aperture position has 
moved significantly at the time-stamp of the candidate 
the candidate is rejected. Of the remaining candidates, 
90% are rejected by visual inspection of the time series. 

Finally we inspect the images of the remaining 23 can- 
didates: they also were all associated with bright moving 
objects overlapping star streaks. No candidate events 
were left in our dataset at the conclusion of this process. 

3.2. Determination of the Stellar Angular Size 

The shape of a lightcurve during an occulta- 
tion event, and hence the detection efficiency, is 
strongly dependent on the angular size of the tar- 
get star (iRoques fc Moncuauetl [2000l iNihei et all 120071 
iBickerton et al.ll2009D . Our fields contain a variety of 
stellar types and a large range of angular sizes (Fig- 
ure |4^). Here we describe the method we use to estimate 
the angular sizes of our target stars in order to account 
for this effect when estimating our detection efficiencies 
(see Section |3"3|) . 

Angular sizes have been related to the position of a 
star in the color - color or color-magnitu de diagrams (e.g., 
Ivan Belld fl999l iNorderen et all l200l. We follow the 
work of INorderen et al.l (2002) and calculate the angu- 
lar siz e of our star targ ets using the 2MASS J and K 
color (jCutri et al. 20031) to invert the set of equations: 

F K = (3.942 ± 0.006) - (0.095 ± 0.007)( J - K) (1) 
Fk = 4.2207 — O.IK — 0.51og# (2) 

where Fk is the surface brightness of a star in K-band, 
which is related to its J — K color, as well as to its 
unreddened apparent K magnitude and angular size 9. 
The relationship between the surface brightness and the 
color of a star (Equation [T]) is calibrated using angu- 
lar s izes measured directly by long baseline interferome- 
trv (|Nordgren et al.|[2002l . 

Not all of our target stars, however, are identified 
2MASS objects, while in the photometry phase we have 
identified all of our targets with USNO-B objects. We 
therefore devised a method that relies on USNO-B R and 
B magnitude to calculate the angular sizes of our targets. 

We first derive the angular size of a subset of targets 
identified with 2MASS objects using the above equations, 
and scale it to obtain the angular size the targets would 
have if their apparent magnitude were R = 12. We 
then considered the USNO-B B — R color for all of these 
targets and calculated a regression on these points. This 
generates a formula that allows us to go from the USNO- 
B color of any of our targets to 2MASS colors and thus 
predict angular sizes according to Equations [T] & [2J for 
an apparent magnitude R = 12. To calculate the true 
angular size we rescale from R = 12 to -Rtjsno 5 - The 
angular sizes of a subset of TAOS targets, rescaled to 
R = 12, is plotted as derived from Equation Q]& [2] versus 
the USNO-B B — R color (Figure |4Jd). Our regression on 
the data is plotted as well (solid line). 

The scatter in the determination of the angular size 
via the method described above is large, as can be seen 
in Figure |4}d. This is due to scatter in the USNO-B color 

5 We do not use our instrumental magnitude for rescaling for 
consistency with what is used in the color determination. 



TABLE 2 
Distribution of synthetic events 



diameter (km) 


implantations 


recoveries 


30.0 


231 


75 


8.0 


385 


84 


3.0 


1078 


73 


2.0 


2003 


89 


1.3 


4393 


73 


1.0 


13255 


66 


0.7 


36222 


40 


0.5 


447764 
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(« 0.3 mag, iMonet et al.l 120031 ) . to the (much smaller) 
scatter in the J and K magnitudes, and to the scatter in 
the empirical determination of the relationship between 
9 and J — K in Equation Q] & [5J We have not used any 
interstellar reddening corrections, and the angular size 
estimation of an unknown reddened star from the near- 
IR relationship would be relatively less affected compared 
to that in visual bands. Reddening is typically small for 
our targets though, since we are only considering objects 
brighter than R ~ 13.5. 

The distribution of angular sizes is well reproduced. 
Figure [4ji shows the distribution of angular sizes for 
a typical TAOS field, calculated via Equations [1] & [21 
and Figure Hh shows the distribution of angular sizes 
in our efficiency simulation obtained via the USNO-B 
color. The distributions do overlap. About 2% of our 
simulated angular sizes fall in the region 9 > 0.15 mas 
and -Rtaos > 11, where there are no observed objects. 
These objects have poor USNO-B color determination. 
Figure HJi shows the region of the 9 — -Rtaos space where 
simulated events are recovered. There are few recoveries 
in the region 9 > 0.15 mas and -Rtaos > d, so these 
stars do not contribute the the expected event rate. 

3.3. Detection Efficiency 

It is necessary to assess the efficiency of our recov- 
ery algorithm in order to derive the number density of 
KBOs from the number of events in our survey. In order 
to measure our recovery efficiency we implant our data 
with synthetic occultations. The data are then repro- 
cessed in the same way we did to search for true events. 
By implanting into the actual lightcurves we do not make 
any assumption regarding the nature of the noise in our 
data. Note that our detection algorithm, described in 
Section [3Tl is not affected by the spectral characteristics 
of the noise, as long as the distribut ion of flux measure - 
ments in a lightcurve is stationary ()Lehner et alj [2010). 
O ur occultation simu lator is based on the work described 
in INihei et all (|2007t l. 

We first generate diffraction lightcurves for KBOs oc- 
culting point sources. We integrate the diffraction pat- 
tern over the disk of our target star. Keeping the stellar 
type fixed, the angular size is modulated by changing 
the apparent magnitude of the star and we can use the 
point source lightcurve to integrate the occultation sig- 
nature over the star disk. A point source lightcurve for a 
D = 3 km KBO at A = 43 AU is shown in FigureEl top, 
and the finite source lightcurve for a magnitude V = 11 
star is shown in the second panel, left. Note the smooth- 
ing of the diffraction features. We modify the lightcurve 
to account for a finite impact parameter b by using the 
finite source lightcurve at impact parameter b = as 






-5 
10 



U 

-5 ^ 
10 - 



(c)i 




-2 



-2 



t (sec) 



t (sec) 



Fig. 6. — Implanted occultations recovered by our pipeline. On the left the lightcurves refer to star targets in a field observed near 
opposition and on the right in a field observed at large angle from opposition, deviation at that epoch from the local mean of the lightcurve 
in units of standard deviation as measured by telescope TAOS D. The value of h\(0) and hg(0), for TAOS A and B, are indicated by 
horizontal arrows labeled A and B, respectively. Each event is described in Tabic \3\ 



TABLE 3 

Parameter of implanted events in Figure [6] 





D (km) 


b (km) 


v le i (km/s) 


SNR 


0* (mas) 


M°) 


M0) 


Mo) 


(a) 


0.7 


0.50 


25.4 


39.8 


0.03 


-4.0 


-6.1 


-5.2 


(b) 


1.0 


0.50 


14.9 


29.4 


0.03 


-5.5 


-6.6 


-7.4 


(c) 


3.0 


1.19 


25.4 


12.8 


0.03 


-5.4 


-7.5 


-8.5 


(<D 


8.0 


2.7 


25.4 


9.4 


0.006 


-6.1 


-12.3 


-7.9 


(e) 


0.7 


0.22 


8.2 


10.2 


0.03 


-5.8 


-3.3 


-6.3 


(f) 


1.0 


0.94 


7.7 


10.8 


0.03 


-10.0 


-8.7 


-8.6 


(g) 


3.0 


1.18 


3.1 


12.9 


0.006 


-2.5 


-4.2 


-8.4 


(h) 


8.0 


2.93 


8.1 


15.7 


0.04 


-11.8 


-13.7 


-11.7 



an input and calculating the intensity of the occultation 
signal at the new distance of each point form the cen- 
ter of the diffraction pattern by interpolating points of 
the finite source lightcurve. A lightcurve for an F0V, 
V=ll star and a 3 km KBO occulting at an impact pa- 
rameter 6 = 2 km is plotted on the right hand side 
of the second panel of Figure [5] Finally, after calcu- 
lating the relative velocity of the KBO as a function of 
distan ce an d angle from opp osition as per iLiang et"al"1 
(120041 ) . and lNihei et~aT1 (|2007f ). we smooth the lightcurve 
to account for finite exposure intervals, and we sample 
the finite exposure lightcurve at the appropriate sam- 
pling rate. In this step we can account for dead-time in 
the sampling interval, which for TAOS is 47.5% at 5 Hz. 
We also allow an offset in time between the center of 
the finite sampled lightcurve and the integration bin. In 
the bottom row of Figure [5] the lightcurves for the event 
in row two are integrated over 105 ms second intervals 
and sampled at 5 Hz, the typical sampling rate of TAOS, 
for an event at opposition and with no time offset (left), 
and for an event at 50° from opposition with an offset 
of 50 ms between the center of the sampling interval and 



the center of the occultation (right). 

In order to sample properly the space of diameters to 
which the survey is sensitive we implant synthetic oc- 
cultations by objects of diameter D = 0.5, 0.7, 1.0, 1.3, 
2.0, 3.0, 8.0, and 30.0 km. For a 30 km diameter KBO 
the event falls in the geometric regime, diffraction effects 
are therefore no longer significant and our efficiency sta- 
bilizes. Because our sensitivity decreases with decreas- 
ing diameter we implant progressively more objects at 
smaller diameters. The number of implantations at each 
size is designed to allow us to obtain a good sampling 
at all sizes. In Table [2] we report the number of objects 
implanted for each size in one of our efficiency runs, and 
the number of recoveries 6 . For objects within the Kuipcr 
Belt (about 30 to 60 AU), the differences induced by dif- 
ferent distances are negligible in the occultation features 
as observed by TAOS. We therefore set the distance to 
A = 43 AU. Every occultation event is implanted at 
a random epoch in the lightcurve set and at a random 
impact parameter between and H/2, where we set H, 

6 Four runs are conducted to improve statistical accuracy in the 
determination of our efficiency. 
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Fig. 7. — Recovery efficiency for 3 km KBOs. Panel (a): efficiency as function of SNR. Efficiency versus crowdedness of the field in panel 
(b), defined as the number of targets brighter than -RtjSNO = 13.5. The efficiency is weighted by the crowdedness. Panel (c): efficiency 
versus magnitude Mtaos ■ I n panel (e) the efficiency as a function of magnitude is weighted by the number of targets at that magnitude. 
Panel (d): efficiency versus relative velocity of the KBO targets. In panel (f) the efficiency versus relative velocity is weighted by the 
relative velocity. All error bars are calculated in a Poissonian fashion from the square root of the number of recoveries. 

a measure of the cross section of the event, to the size 
of the Airy r ing and the proje cted size of the star in 
accordance to iNihei et~all (|2007[ ). 

In order to implant the synthetic occultations into our 
data we modulate the lightcurve by subtracting (adding) 
the amount of flux suppressed (augmented) by the occul- 
tation at each data-point as done in Z08. This approach 
slightly overestimates the noise due to Poisson statistics 
where the flux is suppressed, giving us a conservative 
estimate of our efficiency. We implant exactly one occul- 
tation in each lightcurve in our dataset. 

We then process our implanted lightcurves as we previ- 
ously did to search for events, namely we filter, rank the 
lightcurves, and evaluate the significance £ of each point 
in each lightcurve. We then remove the false positives as 



described in Section I3TT1 Our efficiency decreases rapidly 
with the KBO diameter: from nearly 33% at D = 30 km 
to 2 x 1CT 5 at D = 0.5 km. Note that at D = 30 km we 
ignore diffraction effects and the occultations are mod- 
eled to suppress the flux completely for several consec- 
utive points, depending on the relative velocity, but our 
efficiency is still significantly less than 100%. Some of 
our lightcurves are too noisy to allow detections. 

A set of synthetic events recovered by our pipeline in 
shown in Figure [6l and the parameters of each plotted 
event are given in Table [H Various parameters affect 
our recovery efficiency. Our efficiency for the recovery of 
D = 3.0 km occulting KBOs is plotted in Figure [7] 

The efficiency as a function of SNR is plotted in Fig- 
ure [7K -a few targets at SNR > 100 are left out of 
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Fig. 8. — Effective solid angle for the first 2 years of TAOS data 
(Z08, empty squares) and for the current 3.75 year dataset (solid 
line and filled squares). 

the plot. The behavior of our efficiency as a function 
of crowdedness, where the crowdedness is defined as 
the number of targets in that TAOS field brighter than 
-Rtjsno = 13.5, is plotted in Figure [7p. The efficiency 
is here plotted multiplied by the number of targets in 
the field, to give a better idea of the implication of this 
parameter for the event recovery. The largest number of 
detections is achieved for more crowded fields. 
The efficiency decreases with magnitude (Figure [7b) 
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Fig. 9. — Model-independent upper limits from the TAOS survey 
(solid dots and solid line). Each point represents the upper limit 
to the number of KBO of that size or larger given the TAOS effec- 
tive co verage. Similarly derived upper limits from Bickcrto n et al.l 
112001) at 1 km (BK W), IWang et all lf2009l 1 at 0.5 km (PS) and 
IBianco et all l|200g i at 1 km and 0.7 km (MMT) are plotte d as 
empty circles. The upper limit from Bernstein ct al. ( 2004, di- 
rect survey) is also plotted (HST). The dashed lines represent the 
95% upper and lower limits from S09. The black square (FGS) 
represents the best fit density for the entire S09 survey, with ltr 
error bars. The estimates of the number of objects in the Clas- 
sical Belt (CB), of P lutinos and Scattered Dis k objects (SD) are 
plot ted, as derived by Lcvison & Duncan (1997), Morbidcllij 1119970 
and lVolk fc~Malh otra ( 2003) respectively, assuming each family is 
the unique precursor of JFCs. 



by about a factor of five between magnitude 9 and 13. 
The dominant effect here is the decrease in SNR, though 
a competing effect occurs since lower magnitudes are 
associated with larger angular sizes, and our efficiency 
decreases with increasing angular size (Figure 2J . Fur- 
thermore there are many more dim than bright stars in 
the sky: Figure 03 shows the efficiency as a function of 
-Mtaos multiplied by the number of TAOS targets at 
that magnitude. The highest number of detections hap- 
pen for stars with Mtaos ~ 12.5. 

Our efficiency as a function of the relative velocity of 
the KBO is plotted in FigureEJi. Observing at a pointing 
where the relative velocity of the KBOs is higher boosts 
the event rate of the survey. Our efficiency, however, 
is larger for smaller transiting velocities, particularly for 
small KBOs for which the time-line of the event is shorter 
than one of our data-points at opposition. Ultimately, 
the effective sky coverage of our survey depends linearly 
on both the efficiency and the velocity (see section [3^4)) . 
The efficiency multiplied by the relative velocity u re i is 
plotted against u re i in Figure [7f. Pointing near opposi- 
tion increases the effective coverage of our survey, and 
thus it increases our event rate, for 3 km KBOs. The 
survey strategy can be optimized at different sizes tak- 
ing into account the size dependent efficiency as well as 
the expected KBO size distribution. 

3.4. Effective Sky Coverage and Upper Limits 

We calculate the effective sky coverage of our survey, 
r2e> as: 

1 x , H(D,e r ) ^cl 
\ A 



n e (D) = 



w{D) 



(3) 



where E* is the exposure of the star target (the dura- 
tion of the lightcurve set), A the distance to the oc- 
culter, H the cross section of the event (Section l3.3|) . 
w(D) the weight factor for that diameter, i.e. the frac- 
tion of lightcurves implanted with occultations by KBOs 
of diameter D, and the sum is carried out only over the 
lightcurves where events are recovered. 

The effective coverage of our survey, which takes into 
account our efficiency, is plotted in Figure [8j for both 
the dataset published in Z08 (empty squares) and for the 
current work (filled squares). The solid line is a spline fit 
to the points. 

4. MODEL INDEPENDENT LIMITS ON THE SIZE 
DISTRIBUTION OF KBOS 

We can use fie to calculate model-independent upper 
limits: at each size for which our efficiency calculation 
was conducted we calculate the number density of KBOs 
from the TAOS dataset as a single-point upper limit. 
These limits are shown as filled dots. We can then inter- 
polate these upper limits with a spline fit, obtaining the 
solid line in Figure GO At each diameter D this represents 
the maximum surface density of KBOs of diameter > D. 

Figure [9] shows model-independent upper limits re- 
ported by the occu l tation sur veys of [B ickcrton et all 
(|2008| ) . IBianco et al.l (f200l ) and lWang et al.l (|2009h . The 
improvement is of nearly an order of magnitude at 700 m 
an d over an ord e r of m agnitude at 1 km fro m the results 
of IBianco et al.l (|2009f ) and at 500 m from IWang et al.l 
(2009). As discussed in Sections □ and H the TAOS sen- 
sitivity is limited to relatively bright stars (V < 13.5), 
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Fig. 10. — Expected number of events in 3.75 years of TAOS 
data from the size distribution presented in S09. Slopes of q = 3.9 
(best fit to the HST/FGS data), q = 3.6 and 4.2 (±lcr from the 
best fit) are used. 

and TAOS has a relatively low efficiency at detecting 
D < 3 km KBOs. The large star exposure of TAOS, 
however, provides a dataset that well compensates for 
the lower sensitivity. 

Li mits from the recent HST/FGS occultation sur- 
vey jSchlichting ct al. 2009j, hereinafter S09) are included 
in Figure [§] S09 reported the detection of one pos- 
sible event at a relatively high ecliptic latitude. The 
event is consistent with an occultation by a KBO with 
D = 1.0 km. The dashed lines in Figure [§] represent the 
95% c.l. upper and lower limits from S09. The TAOS 
upper limits are more constraining than the correspond- 
ing limits from S09 for D > 0.6 km. The black square 
(FGS) represents the best fit density for the entire S09 
survey, with la error bars. Note that this is not a model- 
independent estimate, as it is derived assuming a straight 
power law distribution for small KBOs. Note also that 
TAOS has a greater sky coverage at D ~ 1 km, where the 
S09 event was detected, by a factor ~ 6. This detection 
reported by S09 is, however, not statistically inconsistent 
with our upper limit. 

Figure [TTJ1 shows the expected number of events for the 
TAOS survey if one assumes the size distributions de- 
rived in S09. S09 estimated a cumulative size distribu- 
tion T, N (D > 0.5 km) = 2.1 x 10 7 deg -2 from their de- 
tection and non-detections, and fit a power-law to their 
data, anchored at ^n(D > 90 km) = 5.4 deg -2 , to yield 
a differential slope of q = 3.9 ± 0.3. Figure [TOl shows 
the expected number of events for the TAOS survey for 
these distributions. We calculate the number of events 
expected in our survey as: 

Nexp = J Id Qc dD ' 

With no detections our survey can rule out distributions 
that predict -/V cxp > 3 at the 95% confidence level (c.l.). 



The best fit S09 model leads to an expected number of 
2.1 detections in our TAOS data. This is a reasonable 
consistency. We can rule out slopes steeper than q = 4.0 
for this model of the size distribution at 95% c.l. 

Figure |9] also shows the estimates of JFC progenitor 
populations. Our constraints on these families are dis- 
cussed in Section^] 

5. OUTER SOLAR SYSTEM COLLISIONAL MODELS 

The collisional and dynamical evolution of the Solar 
System shaped the size distribution of the Kuiper Belt. 

The belt was originally populated by very small dust 
grains, with small orbital eccentricities (e < 0.01) such 
as tho se we observe in circumst ellar disks around other 
stars (Moro-M artin et "aLI |2007|). Initially thes e small 
objects merge and grow (jKenvon fc Luul[l999a| ): 1 km 
KBOs in the Kuiper belt are thus formed. As their 
gravitational cross section grows larger than their geo- 
metric cross-section, gravitational focusing speeds up the 
growth rate of the largest bodies. This phase is referred 
to as runaway growth, and objects as large as hundreds 
of kilometers can form. One such population, shaped 
primarily through agglomeration processes is predicted 
by theory to have a power law distribut ion in diameter 
dN / dD oc D~ qL with power q^ 4.5 (jKenvon fc Luul 
1999b). Direct observations of large KBOs confirm the 
power law behavior in this regime, the gravitationally- 
dominated region of the size sp ectru m, with a best fit of 
q L = 4.8 (|Fraser fc Kavelaarsl [200l iFuentes fc Holmanl 
120081 and reference therein). The size distribution of 
these large objects, for which gravity dominates the in- 
ternal strength, is remarkably insensitive to parameters 
such as Neptune stirring or the internal tensile strength 
of the KBOs. 

Meanwhile, very large objects in the planetary re- 
gion of the Solar System are also forming into plan- 
ets, t hat are believed to undergo significant migra- 
tions (jTsigani s et al.ll2005l and references therein). The 
orbits of the planetesimals are then stirred up via gravita- 
tional interaction to velocities such that further impacts 
will result in the disruption of the s maller objects: this 
is the catastrophic collisions phase (jDavis fc Farinellal 
fl997t iKmvon fc Luulll999at iMorbidelli et al.ll2008fl : the 
time scale to reach this ph ase is estimated to be be - 
tween 10 Myr and 1 Gyr (jKenvon fc Bromlevl [200i] ). 
For very small objects (probably tens of meters and 
smaller) the collisionally evolved population transitions 
to a regime where the KBO binding energy is dominated 
by internal strength, rather than gravity. Here the col- 
lisional cascade will generate a size distribut ion which 
follows a power law with i ndex q§ = 3.5 (Dohnanvi 
[19691: iKenvon fc Lin J U??99ah . also in a fashion that is 
largely independent on the details of the evolution of 
the protoplanets. Note that the study of collisions be- 
tween icy bodies is still in its infancy, and future work in 
this field will permit assessing the behavior of colliding 
small strength-less or loosely bound particles (jLeinhardtl 
2008). Similarly, future work on coupling collisional 
and dynamical evolut i on co des, recently pioneered by 
ICharnoz fc MorbidelTil (|2007l ). should provide further in- 
sight in the behavior of the size distribution. 

We will refer to the region in between these two regimes 
as the intermediate region. The extent of, and behavior 
of the size distribution in the intermediate region are 
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Fig. 11. — P S05 model and expecte d event rate for TAOS. Left. Triangles are the data from Frascr & Kavelaars (2009). Empty circles are 
the data from Bernstein ct al. (2004). We use these data from direct imaging surveys to set the location of the large end size distribution. 
The model is parametrized with two slopes, q^ = 5 and qi = 3. The positions of the break are plotted at =80, and 20 km: limiting 
values for the PS05 models. Right. Expected yield of events in 3.75 years of TAOS data. The horizontal dashed line represents the highest 
number of events allowed given no detections in the TAOS data (< 3 events at the 95% c.L). Any model above this limit is ruled out by 
our survey. 



instead very sensitive to the formation and evolution pa- 
rameters, and observational information on this region 
can be compared to evolution models. In general, for a 
weaker KBO population the transition to the qs = 3.5 
power law behavior will occur earlier, reducing the size 
of the transition region. A strong KBO population will 
display an extended intermediate region, generally show- 
ing here oscillations around a mean power law of slope 
which also depends on the details of the population and 
its evolution. In this section, we present four models, 
from literature, of the KBO size distribution dN/dD, 
and a simple parametric model. On the basis of these 
models we calculate the number N exp of events expected 
to be detected by the TAOS survey. Since no events were 
found in our survey, any model which predicts N cxp > 3 
is ruled out by TAOS at the 95% c.L 

We caution the reader that the predictions presented 
here depend on both the shape of the size distribution 
and on scaling parameters. All of our models are scaled 
so that the differential size distribution is consistent with 
the direct observations at D = 200 km, where the direct 
surveys are most constraining: 



is constant through the collisional processes. They as- 
sume for most of their model calculation that the inter- 
nal strength of the objects is negligible (gravity domi- 
nated objects). This assumption is motivated by studies 
of comets and asteroids (PS05 and references therein). 
The transition to the fully stren gth-dominated re gime, 
where the size distribution follows lDohnan vi ( 1969) with 
a power law with slope qs — 3.5, occurs at D < 300 m. 
This region of the size spectrum is entirely below the 
sensitivity of TAOS. 

PS05 derive an analytical double power law size distri- 
bution for objects D > 300 m: 



dN/dD 



for D > D b 



dN/dD oc d~ qi for D < D h 



(5) 



dN/dD(D = 200 km) 



0.04 deg" 2 km _1 



(4) 



For each model we will point out the cumulative 
number of KBOs predicted at D — 100 km, I^n(D > 
100 km), with this scaling choice. Strengths and weak- 
nesses of all of these models are also discussed in iFraserl 
(|2009l) . 



5.1. 



Pan & Sari (2005) 

iPan k Saril (|2005l hereinafter PS05) derived a fully an- 
alytical model for the size distribution of KBOs by as- 
suming the population is in a steady state and the mass 



This model is shown in Figure [TTJ left. The slope q has 
value 5l = 5 for large objects and q\ = 3 for objects in 
the intermediate region. PS05 are thus able to calculate 
self-consistently the location of the break in the power 
-Db, which represents the size of the largest KBOs that 
experienced catastrophic collisions, as a function of time. 
The location of the break moves toward larger objects as 
the size distribution evolves. 

Model distributions with break points Db = 80, and 
20 km, limiting values for PS05, and E N (D > 100 km) ~ 
32 deg -2 are shown in Figure [IT] The corresponding 
predicted number of events for the TAOS data analyzed 
in this paper is plotted in the right panel . 

structure of the size distribution in the intermediate 
region is however generally more complicated. Between 

7 Note that a simple extension of the large-size power 
law distribution , whic h represents the naive expectation before 
IBernstein et al] l|2004D appeared, would lead to thousands of de- 
tections. This was the initial design target of the TAOS project. 



13 



R magnitude 
30 35 



M 


OJ 


TZ 




E 




cf 




a' 


in 




1 


•z, 


Q 










CiJj 


c 






D(km) 

Fig. 12. — Our modeling of the KB04 results. The range of 
results of the simulations of KB04 is represented by the shaded 
region (see Figure 10 and 11 in KB04). A few of our models are 
plotted as solid lines. All models are normalized by the slope of 
the large KBO size distribution D~ 4 - 5 . 

Db and the second break, which is the region that TAOS 
can probe, a realistic size distribution is expected to have 
an oscillatory behavior which in the PS05 models pre- 
serves an average slope of q\ = —3 (PS05). These waves 
however, as they are described in PS05, are small in am- 
plitude and affect our prediction to less than 1%. On 
the basis of our no-detection result, break diameters of 
Z?b < 51.3 km are excluded at 95% c.l. This is consistent 
with the data from direct observations and with the au- 
thors' interpretation of the model: in absence of stirring 
by Neptune, the location of the break is consistent with a 
KBO population comprised of objects with little internal 
strength. Note that direct surveys also suggest that for 
one such distribution the location of the break shou ld be 
at a large diameters (see iFraser fc Kavelaara (|2009( ) and 
Figure [TD left). 

5.2. Kenyon & Bromley (2004) 

iKenvon fc Bromley] ()2004 hereinafter KB04) devel- 
oped numerical models of the collisional evolution of 
KBOs, including the effects of the internal strength and 
gravitational binding of KBOs, as well as the initial mass 
in the belt and a model of stirring by Neptune: if in 
its migration Neptune reaches its current location early 
its stirring effect will influence the shaping of the size 
distribution of the Kuiper Belt. The disruption energy, 
defined as the energy necessary to remove 50% of the 
combined mass of the colliding bodies, is modeled as 
Qd = Qb r ^ h + P Qgf^ 9 ; where Qb and /3j, describe the 
internal binding energy and Q g and (3 g the gravitational 
energy, with p the density and r the radius of the object. 
In the KB04 simulations Qb is varied between 10 and 
10 s erg g _1 , p Q g between 10 -4 and 10 4 and (3 g between 
0.5 and 2.0, while /3b is set to 0, as this parameter has 
little effect on the simulation results. 

All KB04 models generally agree in the shape and slope 
of the size distribution for large KBOs (D — 80 km and 
larger) generating a cumulative size distribution which 
follows N(D > 80 km) oc D 35 , equivalent to a power 
law differential size distribution with power q^ ks 4.5, in 
good agreement with data and theoretical predictions for 



the evolution in the gravitationally-dominated regime. 
The models display a variety of behaviors for smaller 
objects. In all simulations a small dip in the 10-40 km 
region is predicted (cf Figures 10 and 11 in KB04). This 
is rou ghly consistent with the results of iBernstein et ail 
(2004): here collisions destroy weak KBOs and models 
with Neptune stirring or weakly bound KBOs produce a 
more significant dip. This feature is followed by an excess 
with respect to the nominal power law for 2-15 km KBOs. 
The amplitude of the excess varies substantially, between 
a factor of a few and a factor of a few tens, depending 
on the internal strength of the objects and on the details 
of the effect of Neptune stirring. The size distribution 
remains sensitive to the details of the models down to 
about D = 50 m, where once again a power law behavior 
begins, with power 3 < q < 5 in the strength-dominated 
regime. 

Figures [12] shows our representation of the models in 
KB04, as they are presented in Figure 10 and 11 of KB04. 
The size distribution is here shown scaled by the slope of 
the large size region (Z? -4 5 ). The shaded region repre- 
sents the range of the KB04 models and the plotted lines 
are some of our models, covering a large region of this 
model space. Focusing on the region near the transition 
between primordial and collisionally evolved population, 
an excess near 5 km and a depletion near 25 km are both 
visible and well represented in our models. 

We model the large size distribution as a power law 
with <7l = 4.5 and the intermediate region with qi = 2. 
We model the excess near D = 5 km with a Gaussian, 
so that: 



dN f , 



dN 
dD 



1 + 7 cx exp- 



(D - p cx f 



(6) 



dD dD \~ ' "~ "' r 2.0 cr e 

We fix the width of the Gaussian excess to a e = 3.5 km, 
we set the location to p cx — 5.5 km or p cx = 1.6 km. 
The intensity of the excess is determined by I ex ; in order 
to fairly represent all results from the KB04 simulations 
we consider models with an excess of / ex = 0, 10 and 
100. Note that as the break diameter moves towards 
large sizes the models naturally simulate the small dip 
near 20 km (Figures IT2"]). Our results are not sensitive 
to the presence of this small depletion. With our scaling 
we obtain a cumulative surface density of KBOs N(D > 
100 km) ~ 25 deg" 2 . 

Figure [T3l shows the models with p, cx = 5.5 km (left) 
and the corresponding expected number of events for the 
TAOS dataset (right) for the limiting break positions 
D = 60 km and D — 10 km for / cx = 0, 10, and 100. 
For a given excess intensity, we can constrain the loca- 
tion of the break: in absence of the excess break diam- 
eters smaller than D^ = 16.5 km are ruled out; break 
diameters D^ < 32.6 km are ruled out for I cx = 10 
and Z?b < 75.3 km are ruled out for an I cx = 100 
excess. When moving the excess towards smaller sizes 
(p-cx = 1.6 km), break locations -Db < 28.3 km are ruled 
out for I cx = 10, and -Db < 63.3 km for I cx — 100. 

Since the KB04 simulations show that models with 
weaker KBOs and Neptune stirring produce a location of 
the break at smaller sizes, our result strongly favors mod- 
els that incorporate the effects of Neptune stirring and 
weaker KBOs, where the bulk strength Qb < 10 3 erg g _1 

5.3. Benavidez & Campo Bagatin (2009) 
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-Dj, = 60 and 10 km the distribution is plotted for a model with no Gaussian excess, and models with Gaussian excess of intensity 7 CX = 10 
and 7 exp = 100, centered at /i ox = 5.5 km. 
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Rec ent simulations by iBenavidez fc Campo Bagatinl 
(|2009t hereinafter BCB09) divide the Kuiper Belt into 
three dynamical families - the CB, the Plutinos (Reso- 
nant Population) and the SD - and follow the collisional 
evolution of each, while taking into account the physics of 
the fragmentation of icy and rocky bodies at the typical 
relative velocities of KBOs. This suite of models ignores 



the effects of Neptune stirring. The models incorpo- 
rate four scaling laws for fragmentation: a simple scaling 
driven by gravitational self-compression, two scaling laws 
which include both self-co mpression and the effe cts of 
strain-rate, as described by iFarinella et al.l (|1982f ). with 
different diameter dependency (D~ - 25 and D~ 5 ) and 
one that follows the modeling of iBenz fc A sphaug (jl999f ) 



15 



R magnitude 
25 30 



R magnitude 



10 - 



O 



2S 
T3 



5Ji 
O 




- 10 



- 1 



0.1 



0.01 



i-3 

> 
o 
w 

CD 
X 
^5 

03 

o 

CD 

PL 

CD 
< 

re 



-= 0.001 



10 



10 2 10 
D(km) 



D(km) 
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for icy bodies. They then vary the material strength of 
the objects between 10 5 and 10 7 erg cm" 3 . 

In BCB09 the size distribution of large objects is set 
to a power law with slope ~ 4.0. This is set as 
an initial condition to their simulations, and since the 
large objects are not undergoing many collisions, this 
shape is preserved. A first break is seen around 100 km. 
The size distribution then departs from power law be- 
havior for all parameter choices, with two to four orders 
of magnitudes fewer KBOs in the ~ 1 km region than the 
nominal power law would predict. The size distribution 
then follows again a power law, with slope qs = 3.5, in 
the strength-dominated regime. The BCB09 models are 
shown in Figure [TJ] left. We simplified the size distri- 
bution behavior in the intermediate region with a single 
power law, though typically the behavior is oscillatory. 
Our scaling leads to N(D > 100 km) ~ 21 deg -2 for the 
BCB09 models. 

These models are all allowed by the TAOS data, 
partly because of the location of the initial break at 
Db = 100 and the slope of the large end size dis- 
tribution. The break location is at the large end of 
what is allowed by direct obser vations (indeed outside 
of the iFraser fc Kavelaarsl ((2008) allowed range of £>b £ 
[50,95 km] with assumption of a 6% albedo). The lo- 
cation of this first break does not evolve in the BCB09 
simulations from what is set as an initial condition. Sim- 
ilarly, the choice of a slope at the large end of the size 
distribution of q^ = 4.0 i s slightly shallower than the cur- 
rent best fit value f rom IFraser fc Kavelaarsl ([2009D and 
iFuentes fc Holmanl (|2008f T As the parameters relative 
to the large size end of the distribution are more firmly 
pinned down by direct observations, these models and 
future occultation data may place stronger constraints 
on the details of the shape of the size distribution below 
the first break and thus the details of the fragmenta- 



tion mechanisms, providing information on the internal 
structure of the KBOs. 

5.4. Fraser (2009) 

iFraserl (|2009l hereinafter F09) considers the collisional 
evolution of the Kuiper Belt size distribution after the 
epoch of accretion. Starting with initial conditions 
that reproduce the observed large end size distribution 
(ql = 4.8), which will not further evolve in the colli- 
sional simulations, and Z?bi = 2 km, the population 
is collisionally evolved over the age of the Solar Sys- 
tem. A depletion, or divot, forms at D ~ 10 x Z?bi, 
or D ~ 20 km. The size of this depletion changes when 
changing parameters relative to the internal strength of 
the KBOs or the impact velocity, as well as choices for 
the initial intermediate slope. An excess is also evident 
at smaller sizes: D ~ 4 km. We model the F09 size 
distributions as they appear in Figure 2 of F09. Our 
size distributions are created as a sequence of 3 slopes: 
<7l = 4.8, 1 < qi < 3, and 2.5 < q s < 4.5, and we model 
the divots and excesses with Gaussians (Figures [T5l left). 
Here we scale the differential size distribution so that 
dN/dD(D = 200 km) - 0.047 deg -1 km -1 , to repro- 
duce the number of objects at D = 50 km 8 . Our survey 
cannot constrain these models: the pronounced divot at 
D 20 km causes a very low event rate in our survey 
(Figure [HI right). 

5.5. Generic 3-Regime Model: Constraints on the 
Intermediate Region of the Size Spectrum 

A simple, generic 3-regime model allows us to describe 
separately the primordial region, the intermediate re- 
gion and the fully collisional, strength-dominated regime. 
Knowing that both the large (gravitationally-dominated) 

8 N(D = 50 km) ~ 15, Fraser, private communication. 
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and the small (strength-dominated) regions of the size 
spectrum are relatively insensitive to the details of the in- 
ternal structure and evolution of the Kuiper Belt, it is the 
intermediate region that contains the most information 
about the physical details of the KBOs. We grossly sim- 
plify the expected structure in the intermediate regime 
and describe it with a power law. The parameters of one 
such models would then be three slopes q^, qi, and qs , 
and two break locations Dbi and -Db2- 

We set the slope of the large end of the size spec- 
trum and the location of the first break to the best 
fit to the data from direct observations: q^ = 4.8 and 
Dbi = 75 km. We model the intermediate region as a 
plateau, or with a shallow slope: q\ £ [0, 3], and the small 
size end of the spectrum as a power law with slope 3.5 
for the strength-dominated, collisionally evolved popula- 
tion. Figure [iUl shows our 3-regime model (left) and the 
limit we can set to the intermediate slope, qi, and sec- 
ond break location, Db2, phase space (right). Any pair 
of values q\ and Db2 that fall above the solid line are 
excluded. 

6. THE JUPITER FAMILY COMETS PROGENITOR 
POPULATION 

The Jupiter Family Comets (JFCs) are believed to 
originate in the Outer Solar System. In this scenario 
the giant planets generate gravitational perturbations 
that affect the orbits of the Outer Solar System bod- 
ies, injecting them into the planetary region, where they 
are captured by Jupiter. The orbital inclination of the 
JFCs suggests that their precursor population has a 
disk-like distribution, favoring t hus the Kuiper Belt over 
the Oort Cloud as a reservoir ()Volk fc Malhotral 120081 
and references therein). The Classical Belt (CB), the 
Plutinos and the Scattered Disk (SD) have been con- 
sidered as precursors in various studies. The dynam- 
ical characteristics of each population determines the 
efficiency of the injection process and the number of 
objects in each progenitor family can thus be derived 
on the basis of the den sity of JFCs, which is observa- 
tionally constrained fsee iTancredi et al"1l2006L and refer- 
ences therein) . Furthermore the size distribution of JFCs 
should re flect the size distribut ion o f the progenitor pop - 
ulation. iBernstein et~aTI (|2004D and lBianco et al.l (|2009f) 
have pointed out that a better determination of the size 
distribution of the Kuiper Belt would help understand 
the origin of the JFCs. 

We assume the JFC precursors are in the size range 
1 — 10 k m as this is obser ved to be the typical size 
of JFCs (|Lowrv et al.l 120081) . We considered the esti- 
mates on the KBO populations (CB and Plutinos) and 
SD derived from dynamical simulations under the as- 
sumption that each population is the unique progeni- 
tor of the JFCs, and we compared them to the upper 
limits deriv ed from our survey (Figur e E} . We use the 
estimate of iLevison fc Duncan! ()1997f ) for a population 
of cometar y prec ursors entirely in the CB, and that of 
IMorbidellil (|1997l ) for Plutino progenitors. These are 
converted into a surface density by assuming for each 
populati on a projected sky a rea of 10 4 deg 2 , as was 
do ne by [Be rnstein et al. (20Q3). We consider the results 
of lVolk fc Malhotral (f200l for a progenitor population 
in the SD . We calculate the minimum surface density 
of SD objects expected in the region of sky typically ob- 



served by the TAOS survey. For this we use information 
on the fraction of time the objects spend between 30 and 
50 AU and within 3° of the ecliptic plane as provided by 
IVolk fc Malhot ra (2008). These estimates on the number 
of objects are shown in Figure [9] as horizontal lines. 

Our results rule out a precursor family composed 
uniquely of CB objects with D > 3 km. Occultation 
surveys are the only surveys that at present can probe 
this region of the size spectrum, and our preliminary re- 
sult shows that future occultation surveys will be able to 
derive useful constraints on the origin of the JFCs. 

7. CONCLUSIONS 

We presented an analysis of 3.75 years of TAOS data, 
comprising 5 x 10 5 star-hours of lightcurves sampled at 
4 or 5 Hz observed with three telescopes simultaneously. 
We searched for occultations of our target stars by KBOs 
in order to constrain the size distribution of KBOs, par- 
ticularly in the 500 m to 10 km region, which is currently 
out of reach of direct observation surveys. More than 
90% of the TAOS data is collected within 5° of ecliptic 
latitude in order to maximize the occultation rate by ob- 
jects in the Kuiper Belt. Occultations near opposition 
lead to a higher event rate for TAOS, even after taking 
into account the increased recovery efficiency for small 
objects where the angle from opposition is larger and 
the relative velocity of the KBOs is lower. 

We found no occultation events in our data. This al- 
lowed us to set upper limits to the number density of 
KBOs that are stringent enough to be compared use- 
fully with models for the formation and evolution of 
the Kuiper Belt. We considered four theoretical mod- 
els, PS05, KB04, BCB09 and F09, all of which describe 
the present size distribution of the Kuiper Belt, and we 
set constraints on these models. This is the first detailed 
comparison of occultation data with specific model re- 
sults. 

Our result, particularly when compared with PS05 and 
KB04, suggests that the Kuiper Belt is populated by 
fragile bodies, and that the effect of the migration of 
Neptune played an important role in its formation. None 
of the BCB09 or F09 models can be ruled out. 

Using a generic model, where the size distribution is 
described by three consecutive power laws, and fixing 
the slope on the large end size and the location of the 
first break to the best fit from direct observations (^l = 
4.8 km, L>bi = 75 km), and the slope of the small side 
to q s = 4.0 we can constrain the intermediate slope q\ 
and the location of the second break £>b2, as shown in 
Figure 1161 As direct surveys are currently not sensitive 
to KBOs smaller than D ~ 28 km occultation surveys 
provide the only probe of this region of the size spectrum, 
and the large TAOS dataset allowed us to set the first 
constraints to the location of the second break. 

We also considered the Jupiter Family Comets. As- 
suming the JFCs are injected into their present or- 
bit from one of the Kuiper Belt populations, Classical 
Belt, Plutinos, or from the Scattered Disk, we compared 
the upper limit derived from our survey to the esti- 
mates of the n umber of objects de r ived u sing the num- 
ber of JFCs by ILevison fc Duncanl (|1997| ) for a popula- 
tion of cometary precursors entirely i n the CB, that of 
IMorbidellil (|1997l ) for Plutinos and of IVolk fc Malhotral 
( 2008f ) for a progenitor population in the SD. We can 



17 



24 



R magnitude 
26 28 30 



32 



o 



2S 
T3 



o 



4 - 



3 - 



2 - 



- 




- 



100 70 



D(km) D b2 (kni) 

Fig. 16. — Three-slope model described in Section 15.51 (left) for the minimum second break location Z?t,2 allowed by slopes qj 
3. Slope— break location phase space for the three-slopes model (right). The shaded region is excluded by our data. 



and 



rule out the a unique precursor family composed of CB 
objects D > 3 km. This preliminary result confirms that 
occultation surveys can help understanding the origin of 
JFCs. 

A recent analysis of the HST guiding data (S09) re- 
ports the detection of an occultation by a D ~ 1 km 
KBO. This object is detected at high ecliptic latiutde 
(b ~ 14°). The surface density derived in S09, assuming 
a straight power law size distribution for small KBOs, 
is within the upper limit set by TAOS, although a la 
increment over the best fit surface density is ruled out to 
better than 95% c.l. 

TAOS has operated for over four years observing con- 
tinuously with three, and now four, 50 cm aperture tele- 
scopes. TAOS is only marginally sensitive to sub-km 
KBO occultations, but we were able to prove that the 
low sensitivity at sub-km sizes is more than compen- 
sated for by the vast exposure of which TAOS is capable. 
Improvements in our selection criteria and the analysis 
of four telescope data should increase our sensitivity to 



smaller objects. 
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